In [54]:
import numpy as np
from math import exp, pow, sqrt
from numpy.linalg import inv
from functools import reduce
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from matplotlib.colors import ListedColormap

Interpolation with RBF

$$ f( x) =\sum ^{P}_{p=1} a_{p} .R_{p} +b $$$$ R_{p} = e^{-\frac{1}{2\sigma ^{2}} .\parallel ( X_{i}) -( X_{p}) \parallel ^{2}} $$$$ \sigma =\frac{P_{max} -P_{min}}{\sqrt{2P}} $$$$ \sigma =\frac{4-2}{\sqrt{2.2}} \ $$$$ \sigma ^{2} =1 $$$$ C_{1}=2 $$$$ C_{2}=4 $$
$$\displaystyle \frac{1}{[ R]} =\left([ R]^{t} .[ R]\right)^{-1} .[ R]^{t}$$
$$ \displaystyle \begin{bmatrix} a \end{bmatrix} =\frac{1}{[ R]} \ \begin{bmatrix} A \end{bmatrix} $$

In [55]:
def rbf(inp, out, center):
    def euclidean_norm(x1, x2):
        return sqrt(((x1 - x2)**2).sum(axis=0))
    
    def gaussian (x, c):
        return exp(+1 * pow(euclidean_norm(x, c), 2))
    
    R = np.ones((len(inp), (len(center) + 1)))

    for i, iv in enumerate(inp):
        for j, jv in enumerate(center):
            R[i, j] = (gaussian(inp[i], center[j]))
            
    Rt = R.transpose()
    RtR = Rt.dot(R)
    iRtR = inv(RtR)
    oneR = iRtR.dot(Rt)
    a = oneR.dot(out)
    
    def rbf_interpolation(x):
        sum = np.ones(len(center) + 1)

        for i, iv in enumerate(center):
            sum[i] = gaussian(x, iv)

        y = a * sum
        return reduce((lambda x, y: x + y), y)
        
    return rbf_interpolation


In [35]:
inp = np.array([2, 3, 4])
out = np.array([3, 6, 5])
center = np.array([2, 4])

rbf_instance = rbf(inp, out, center)

In [36]:
input_test = input_test = np.linspace(0,10,100)
output_test = list(map(rbf_instance, input_test))

In [37]:
plt.plot(input_test, output_test)
plt.plot(inp, out, 'ro')
plt.ylabel('expected vs predicted')
plt.savefig("rbf1.svg")
plt.show()




In [63]:
inp = np.array([2, 3, 4, 5])
out = np.array([3, 1, 5, -2])
center = np.array([2, 3, 4])

rbf_instance = rbf(inp, out, center)

In [64]:
input_test = np.linspace(-5,10,100)
output_test = list(map(rbf_instance, input_test))

In [65]:
# plt.plot(input_test, output_test)
plt.plot(inp, out, 'ro')
plt.ylabel('expected vs predicted')
plt.savefig("interpolate1.svg")
plt.show()




In [60]:
inp = np.array([2, 3, 4, 5])
out = np.array([3, 1, 5, -2])
center = np.array([2, 3, 4])


rbf_instance = rbf(inp, out, center)

In [58]:
input_test = input_test = np.linspace(-5,15,100)
output_test = list(map(rbf_instance, input_test))

In [59]:
plt.plot(input_test, output_test)
plt.plot(inp, out, 'ro')
plt.ylabel('expected vs predicted')
plt.savefig("rbf3.svg")
plt.show()



XOR input


In [44]:
inp = np.array([np.array([1,1]), np.array([0,1]), np.array([0,0]), np.array([1,0])])
out = np.array([              0,               1,               0,             1])
center = np.array([ np.array([1,1]), np.array([0,0])])

rbf_instance = rbf(inp, out, center)

In [45]:
inp_test = np.array([np.array([1,1]), 
                     np.array([0,1]), 
                     np.array([0,0]), 
                     np.array([1,0])])
output = map(rbf_instance, inp_test)

In [46]:
def colorize(output):
    c = [None]* len(output)
    for i, iv in enumerate(output):
        if (output[i] > 0):
            c[i] = 'blue'
        else:
            c[i] = 'red'
    return c

In [53]:
inp_x = [1, 0, 0, 1]
inp_y = [1, 1, 0, 0]

c = colorize(output)

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

xx, yy = np.meshgrid(np.arange(0, 1, 0.02), np.arange(0, 1, 0.02))

cm_bright = ListedColormap(['#FF0000', '#0000FF'])

ax.scatter(xx[:, 0], yy[:, 1], output, cmap=cm_bright, depthshade=False)
plt.savefig("rbf_xor.svg")
plt.show()


---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-53-f884d345d84f> in <module>()
     12 cm_bright = ListedColormap(['#FF0000', '#0000FF'])
     13 
---> 14 ax.scatter(xx[:, 0], yy[:, 1], output, cmap=cm_bright, depthshade=False)
     15 plt.savefig("rbf_xor.svg")
     16 plt.show()

/usr/local/Cellar/python@2/2.7.15_1/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/mpl_toolkits/mplot3d/axes3d.pyc in scatter(self, xs, ys, zs, zdir, s, c, depthshade, *args, **kwargs)
   2354 
   2355         xs, ys, zs = np.broadcast_arrays(
-> 2356             *[np.ravel(np.ma.filled(t, np.nan)) for t in [xs, ys, zs]])
   2357         s = np.ma.ravel(s)  # This doesn't have to match x, y in size.
   2358 

/usr/local/lib/python2.7/site-packages/numpy/lib/stride_tricks.pyc in broadcast_arrays(*args, **kwargs)
    247     args = [np.array(_m, copy=False, subok=subok) for _m in args]
    248 
--> 249     shape = _broadcast_shape(*args)
    250 
    251     if all(array.shape == shape for array in args):

/usr/local/lib/python2.7/site-packages/numpy/lib/stride_tricks.pyc in _broadcast_shape(*args)
    182     # use the old-iterator because np.nditer does not handle size 0 arrays
    183     # consistently
--> 184     b = np.broadcast(*args[:32])
    185     # unfortunately, it cannot handle 32 or more arguments directly
    186     for pos in range(32, len(args), 31):

ValueError: shape mismatch: objects cannot be broadcast to a single shape

In [ ]: